Generation of short and long range temporal correlated noises 
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We present the implementation of an algorithm to gener- 
ate Gaussian random noises with prescribed time correlations 
that can be either long or short ranged. Examples of Langevin 
dynamics with short and long range noises are presented and 
discussed. 



I. INTRODUCTION 

In stochastic process simulations, the generation of 
noise with some specific statistical properties has been 
of great importance. In most cases, the required noise 
is Gaussian and white (delta correlated); in many other 
cases, as in noises of real systems, a specific correlation is 
needed and an appropriate noise generation has to be im- 
plemented which can be carried out in standard computer 
facilities. Despite the extensive work in developing al- 
gorithms to generate noises with particular correlations, 
there is still a lag concerning those noises with any given 
temporal, spatial or spatiotemporal correlation. 

Some algorithms have been proposed in the last few 
years for noises which have been proved to obey a lin- 
ear Langevin equation with a linear, Gaussian, white, 
noise term. A formal integration of this linear equation 
is enough to generate a random process with a very par- 
ticular time correlation. Examples of this approach are 
the so called Ornstein-Uhlenbeck (O-U) and Wiener (W) 
processes. The physical meaning of the O-U process is 
the velocity of a Brownian particle under the influence of 
friction and immersed in a heath bath. W-process is the 
archetype for the dynamical behavior of change in posi- 
tion of a free Brownian particle in the high friction limit. 
The correlations of these two processes are well known 
and can be obtained analytically |Q and numerically j^] 
by standard methods of stochastic processes. 

Nevertheless, quite often in numerical stochastic sim- 
ulations, we are faced with processes characterized by 
a Gaussian noise with a specific correlation and an un- 
known Langevin-like equation dynamics. Therefore, a 
more general algorithm that only depends on the knowl- 
edge of the temporal correlation is necessary . 

The purpose of this work is to present the implemen- 
tation of an algorithm which allows to simulate Gaussian 



noises with almost any given temporal correlation func- 
tion. The only requirement for the algorithm to work 
is that the Fourier transform of the temporal correlation 
function should be known Q . In cases where the Fourier 
transform of the required temporal correlation function 
can not be properly defined, the algorithm can be used 
with a suitable cutoff. 

In Section II, algorithms based in Langevin equations 
are revised and the new algorithm of noise generation 
is presented in practical terms. The algorithm has the 
property of simulating different correlation regimes, from 
short to long ranged. Section III considers different ap- 
plications for which the dynamics of the system are very 
sensitive to the noise correlation. A summary of results 
and some conclusions are presented in the last section. 

II. ALGORITHMS TO GENERATE GAUSSIAN 
NOISES 

A. The Linear Langevin Equation Method 

For the sake of comparison with the method we want 
to present here, let us revise briefly the main points of 
the generation of a Gaussian noise which follows a linear 
Langevin equation in terms of a Gaussian white noise. 
The stochastic discretized trajectories are generated by 
formally integrating the corresponding Langevin equa- 
tion and using a set of Gaussian random numbers for 
a given realization. The statistical properties are calcu- 
lated from many realizations of these trajectories ("an 
ensemble"). In the case of Gaussian noise, only two mo- 
ments are necessary to characterize the statistical prop- 
erties of the random process. Due to the linear character 
of the Langevin equation the Gaussian property of the 
white noise is transmitted to the generated noise. 

The simplest case is the W-process, which follows the 
Langevin equation 

v = ew , (i) 

where is a Gaussian white noise of zero mean and 
delta correlation: 

mm) = 2eS(t-t') . (2) 
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Here, e is the noise intensity. 

The statistical properties are easily evaluated fw and 
a trajectory for this type of noise can be obtained by 
formally integrating Eq.(Q) over time: 

pt+At 

V (t + At) = n(t) + J Z(t')dt' = m + x(t) , 

(3) 

where the noisy term X(t), is constructed from 



X(t) = v / 2eAi a(t) 



(4) 



Here, a(t) is a set of Gaussian independent random num- 
bers of zero mean and variance equal to unity, obtained 
from any reliable Gaussian random generator . 

As a second example, we consider the O-U process. 
This process contains a temporal memory and obeys the 
Langevin equation 
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v + m 

T T 
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where r is the characteristic time memory. As in the 
former formal integration gives: 

r)(t + At) = T](t)e-^ H / e- 1 -^- dt' . 

(6) 

Studying the statistical properties of the noisy term, the 
algorithm reads M , 



V (t + At) = r,(t)e-- + ^-(l-e-f)a(t). (7) 

In this way, a time stochastic trajectory is generated 
step by step as in the former example. Moreover, in the 
algorithm we introduce below, the whole trajectory is 
constructed in a single calculation step. 



B. Spectral method 

As wc noted in the introduction, this method starts 
from the knowledge of the time correlation function. 
Probably the basis of this algorithm has been rediscov- 
ered and implemented several times, but now, due to the 
actual high speed and large storage space in computers, 
the algorithm can be implemented rather easily. Here 
we follow the main ideas of Ref. HQ and first introduce 
some steps that simplify the calculations. 

We want to generate a Gaussian, random correlated, 
noise, rj(t), whose correlation function j(t), is defined by: 



( V (t)v(t')) = 7(1* -f I) 



(8) 
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is known (to some extent). In the w- Fourier space, this 
correlation reads: 



{r)(u)r}(u')) = 2ttS(uj + w') j(uj) 



(10) 



where rj(uj) is the Fourier transform of r)(t). 

With this initial information in mind, the algorithm 
could be summarized as follows. First we discretize the 
time in N = 2" intervals of mesh size At, and note that 
this time interval has be much smaller than any other 
characteristic time of the system for our method to work. 
Every one of these intervals will be denoted by a Roman- 
index in real space (time) and by a Greek-index in Fourier 
space (frequency). 

In the discrete Fourier space, the noise has a correla- 
tion given by: 



Now, f?(w M ) can be constructed from 



(11) 



/i = 0...A, r](u> ) = ii(uj n ) , 



W/1 = 2*£ (12) 
p A At ' v ' 



where a(w^) = a M are Gaussian random numbers with 
zero mean and correlation, 



= S. 



(13) 



This type of delta-anticorrelated noise can be generated 
rather easily if the symmetry properties of real periodic 
series in the Fourier space {ai) are used Q. Thus, we 
avoid the extra work involved in Fourier transforming 
real random numbers Q . Consider for example a system 
of size N, the Fourier components of a periodic series are 
then related by: 



Otfji+pN 



a,, 



a_ 
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where p is an integer number. Since a^ois real and 
Im(a fl> o) = — 7m(a M< o), as can be seen in figure [l], 
the requested anticorrelated complex random numbers 
(|J|) can be constructed as 



ib,, 



b 







(15) 



where and 6 M represent Gaussian random numbers 
with zero mean and a variance of one half: 



1 



(16) 



and its Fourier transform, 
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FIG. 1. Diagram representing the construction of a discrete 
delta anticorrelated noise in Fourier space. Indexes a; and hi 
represents Gaussian random numbers as explained in the text. 

The discrete inverse transform of ?/(wp) is then nu- 
merically calculated by a Fast Fourier Transform algo- 
rithm |6|. The result is a string of N numbers, r)(U), 
which by construction, have the proposed time correla- 
tion (|8|). However, due to the symmetries of the Fourier 
Transform, only N/2 of these values are actually inde- 
pendent and the remaining numbers are periodically cor- 
related with them. 

In order to check the suitability of the procedure, the 
time correlation of Equation (||) is numerically evaluated 
by an independent (non Fourier based) method, namely, 



7(*i) 



N„ 
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(17) 



Here, N max is a number smaller than N /2 (in the exam- 
ples that follow, it is taken equal to N/i). 

We now present several examples where this approach 
is implemented. Since the number of applications making 
use of temporal random noises is quite large we have cho- 
sen applications with rather different correlation proper- 
ties in order to illustrate the power of the method. The 
same procedure has already been used for spatial corre- 
lated noises @,|,[ii). 



C. Short range correlated noises 

a.l. A Gaussian Correlation. 

Let us consider a noise with a Gaussian correlation 
function defined as 



(v(t)v(t')) = i(\t-t'\) = 



2c 



Ty/2^ 



(18) 



where e and r are the noise intensity and correlation time, 
respectively. The correlation is normalized in such a way 
that 



j(t) dt 



(19) 



Setting t — > 0, the white noise limit is recovered. The 
Fourier transform of this Gaussian correlation(|l8|) is 
given by: 



7(0;) = 2ee 4 



(20) 



According to the prescription ([L2|), we now generate 
the discrete field ?y(w M ) as: 



77(0^) = NAt 2ee 



,±p(.cos(.2wu/N)-l) 



1/2 



(21) 



In Fig an explicit comparison between the numeri- 
cal results and the expected theoretical prediction is pre- 
sented. 




FIG. 2. Temporal correlation for two short range cor- 
related noises. a) Gaussian correlation b) 
Orstein-Uhlenbeck process (p^j. Common parameters: 
N = 2 17 , At = 0.01, e = 20, r = 10. Full lines are the simula- 
tion data and dashed lines are the corresponding theoretically 
expected results. 

a.2. The Ornstein-Uhlenbeck process. 

The Ornstein-Uhlenbeck process simulates the behav- 
ior of the velocity of a Brownian particle under friction 
and immersed in a thermal bath. Quite often, it is used 
to represent a real noise with memory, whose intensity is 
e and whose correlation time (or memory intensity) is r. 
This is a well known Gaussian, Markovian, and station- 
ary noise, which obeys a linear Langevin equation (|^). 
Its well known correlation is given by: 



(v(tW)) = i(\t-t'\) 



e it-t'i 
—e 

T 



(22) 



with the same normalization used for the previous exam- 
ple. Since this noise obeys a linear Langevin equation, 
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which was not the case for the first example, a realization 
of this noise could be simulated using the formal solution 
of the stochastic differential equation given byEq.(@). In- 
stead, we present in this subsection the results obtained 
following our spectral method. The Fourier transform of 
this particular correlation function is: 



7M 



2e 



1 + {TLjf 



(23) 



According to the prescription (|l2|), we generate a dis- 
crete field rj{uj^) as: 



NAt 2e 



1 + (^«n(7rp/A0) 2 



1/2 



(24) 



In Fig. ||, we compare the numerical and theoretical 
results. The discretization of the w-variable using ei- 
ther the function cos(^ i ) as in Eq.(pT|) or the function 
sin(^) as in Eq. (|24|) does not make any difference. 

We have seen that in the short range correlated noises 
considered here, the agreement between the statisti- 
cal properties of the noise generated from the spectral 
method and the statistical properties required from the 
noise are quite good. 



D. Long range correlated noises 

A more dificult noise, in the sense that it can not be 
obtained from any known linear Langevin equation, is 
the one characterized by a power-law decaying correlation 
function, 



7(1* -A) 



\t-tr 



<p< 1. 



(25) 



This correlation is not well defined in Fourier space (it 
has a singularity at \t — t'\ = 0). Therefore, in order to 
implement the spectral method, we start with a guess for 
7(W/i), and then we look at the dynamics of j(t) in real 
space. The Fourier values of the noise are now discretely 
generated according to the expression, 



NAt 



•wo 



|(/3-l)/2 



(26) 



where ujo is a predefined cut-off. In Ref. m modified 
Bessel functions are used for the correlation in Fourier 
space but they have the same long range decay. In our 
case, we assume the following form for the correlation 
function in real space, 



7(t) 



Ac 



< < 1. 



(27) 



where t = At/it, (3 is the parameter describing the power 
law decay, e is the noise intensity, and A is a parameter 
to be fitted from a correlation average. 



Fig. ||, presents two examples of the power-law Gaus- 
sian noise with j3 = 1/3 and /3 = 2/3. It can be seen that 
the numerical results are fitted very well by the expected 
power laws. 




FIG. 3. Temporal correlation for two long range corre- 
lated noises, with correlation defined by equation ^2o\). a) 
P = 1/3, e = 20,A = 0.7095. b) = 2/3, e = 20, A = 0.4483. 
Common parameters: N — 2 17 , At — 0.01. Full lines are the 
simulation data and dashed lines are the corresponding theo- 
retical fitting from equation §P\). 

The existence of power law noises has been predicted 
and discussed for quite some time in the literature. Many 
examples arise from theoretical as well as from exper- 
imental studies. A recent reported case is the exami- 
nation of the temperature fluctuations in climatological 
data |l2|,[l3| , where a power law functional form is found 
for the correlation of these fluctuations. In this partic- 
ular case, the analysis is carried over a temperature au- 
tocorrelation function as defined by Eq. ( 17), and it is 
found that the distribution of temperature fluctuations is 
well described by a Gaussian distribution with the long- 
ranged decaying auto-correlation function: C(s) ~ s -7 . 
The exponent found for their data is around 2/3. Noises 
with these properties can be generated using our pro- 
cedure and can be implemented in modelling processes 
similar to the one reported in this climatological anylisis. 
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FIG. 4. Section of the time series realization for a power- 
law correlation with (3 = 2/3 

Fig. ^ shows a realization of a particular noise trajec- 
tory generated with our method for a power law noise 
with (3 = 2/3. The variations of the data generated with 
the algorithm look very similar to figures l.b, 6. a and 6.b 
in ref. |l^] . The persistence of the series is evident from 
our simulation where the stochastic trajectory appears 
as "packets". 



range noises. The explicit form of the time dependence of 
this quantity appears in Fig. ||, where we also show the 
variance of the three different examples in which our al- 
gorithm is applied. Two super-diffusive cases are clearly 
seen, with dispersion exponents ~ 1.75 and 1.25, corre- 
sponding to noise decay power laws, (3 — 0.25 or 0.75, 
respectively. A short range noise (exponential) is also in- 
cluded for comparison. We see that in this last case the 
behavior is ballistic (deterministic) ~ t 2 for t < r, but 
diffusive for t 3> r as predicted by Eq. (g9). 
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III. APPLICATIONS 

We proceed now to present and discuss two differ- 
ent non-trivial examples in non-equilibrium statistical 
physics where our algorithm can be applied: (i) the dis- 
persion process of a Brownian particle and (ii) the decay 
from an unstable state, both cases under the influence of 
long range noises generated through the spectral method. 
In the first case the noise is additive whereas in the sec- 
ond case it appears multiplicatively. 

A. Superdiffusive motion 

In this section we want to focus our attention on the 
random motion of a Brownian particle which obeys the 
Langevin equation 

§ = *> • < 28 > 

where r\ is a Gaussian noise with a given correlation func- 
tion. We consider three different cases, a short range and 
two long range noises. 

The solution for the relative dispersion is well known: 

(Sx{t) 2 ) = 2 f f y(s) ds dt'. (29) 
Jo Jo 

For large times, this expression is either ~ t (diffusion) for 
short range noises, or ~ t 2 ~P (super-diffusion) for long 



FIG. 5. Dispersion as a function of time for three different 
realizations of the noise. The two upper curves correspond to 
a power law noise with two different exponents, (3 — 0.25 and 
(3 — 0.75. The bottom curve corresponds to an exponential 
noise with r = 10.0. Straight lines are plotted as an eye guide. 

B. Decay of an unstable state 

The decay of an unstable state is one of many interest- 
ing problems which appear in noncquilibrium phenomena 
and nonlinear relaxation process. The switch-on process 
of a dye laser has been a prototype of a noncquilibrium 
situation in which the influence of several sources of noise 
have been tested. The dynamics of this system can dis- 
criminate the effects of both additive or multiplicative 
noises JT^-p^l . It has been established that additive noise 
is responsible for the short time dynamics but multiplica- 
tive noise effects appear in the medium and long term dy- 
namics. In previous studies, white or short range, colored 
noises have been used, but long range noises have never 
been considered. Here we will present several numeri- 
cal simulation results for the decay of an unstable state 
under the influence of long range multiplicative noises. 

In order to simplify the analysis, let us illustrate the 
situation with the following Langevin equation, 

^ = ax - bx 3 + xritt) , (30) 
dt 

where r/ is a Gaussian noise whose correlation function 
is to be specified. For a and b positive parameters, the 
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initial value x(0) = is an unstable state which does need 
small perturbation to start relaxing towards its steady 
state. To trigger this process we need either the additive 
noise of Eq. (|30| ) or an initial distribution for x(0). We 
have chosen this last assumption and the initial values of 
x(0) are Gaussian distributed with statistical moments 
(x(0)) = and (x(0) 2 ) = a 2 . Due to the symmetries of 
the problem, the mean value is (x(t)) = 0, and hence we 
look at the dynamical evolution of the second moment. 

To get a precise idea of the short time behavior, it 
is enough to study Eq. ( |30| ) in a linear approximation. 
Formal integration of this equation gives, 

x(t) = x(0) exp (at + r]{t') dt'^j . (31) 

Using Gaussian properties through the calculation, we 
obtain for the second moment 

(x 2 (t)} = a 2 exp{2at + 4fi(t)) , (32) 

where f2(£) is defined by: 

ft(t) = [ ( 7(s) ds d£ . (33) 
Jo Jo 

For the deterministic case, e — 0, and the multiplica- 
tive white noise case, we have: 

(x 2 (t)) = a 2 exp(2at) , (34) 

(x 2 (t)) = a 2 exp (2at + Aet) , (35) 
respectively. 

As in the superdiffusion case, we get for long range 
noises a power time dependence ~ t 2 ~P in the exponen- 
tial; we therefore expect a larger rate for the relaxation 
of the initial state. 




In Fig. |6| we see how the decay of the unstable initial 
state is influenced in the different cases discussed. For the 
values of the parameters we have used, Eqs.(|34|) and fl35| ) 
predict that the multiplicative white noise should decay 
six times faster than that the deterministic case. For long 
range correlated noises we expect that the smaller the 
exponent (3 is, the faster it will decay. Within statistical 
errors, these points are actually seen in the simulation 
data presented in the figure . 

It is interesting to note that for intermediate and long 
times no analytical results can be obtained. However, 
figure ^ is telling us that the final steady state also de- 
pends on the correlation of the noise, taking larger values 
for smaller /3's. This is not an easy problem and it would 
need further theoretical analysis. 

IV. SUMMARY AND CONCLUSIONS 

We have presented a method to generate Gaussian 
noises with a prescribed time correlation function which 
does not depend on any type of dynamics. The algo- 
rithm herein incorporates whithin the same framework 
the generation of long and short ranged noises. In par- 
ticular, the need for a good algorithm is manifest when 
long range noises representing physical dynamical process 
with unknown or very complicated dynamical equations 
are required. In this sense, this method can numerically 
simulate the influence of long and short ranged noises in 
physical systems in a very reliable and controlled way. 
With the prescription given, many problems with long 
range realizations can be simplified and the possibility 
of some kind of analytical approach to specific problems 
increases (such is the case for the examples we have pre- 
sented) . 

The generalization of the algorithm to spaciotemporal 
noises is straight forward for any number of dimensions. 
Since by construction, the noise generation is finite, the 
upper time limit of the simulation has to be selected in 
advance, but it does not constitute a serious limitation 
for the algorithm. 
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